Elliptic functions and efficient control of Ising spin chains with unequal couplings 
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In this article, we study optimal control of dynamics in a linear chain of three spin 1/2, weakly 
coupled with unequal Ising couplings. We address the problem of time-optimal synthesis of multiple 
spin quantum coherences. We derive time-optimal pulse sequence for creating a desired spin order 
by computing geodesies on a sphere under a special metric. The solution to the geodesic equation 
is related to the nonlinear oscillator equation and the minimum time to create multiple spin order 
can be expressed in terms of an elliptic integral. These techniques are used for efficient creation of 
multiple spin coherences in Ising spin-chains with unequal couplings. 

PACS numbers: 03.67.-a, 82.56.-b 



I. INTRODUCTION 

In the absence of relaxation, experiments in coherent 
spectroscopy and quantum information processing con- 
sist of a sequence of unitary transformations on the quan- 
tum system of interest. Pulse sequences that create a 
desired unitary transformation in the minimum possible 
time are of particular interest in experimental realiza- 
tions as they can reduce losses due to relaxation. This 
poses the problem of time optimal control of quantum 
systems, which is of both theoretical and practical in- 
terest in the broad area of coherent control of quantum 
systems. 

Besides application in spectroscopy, synthesizing uni- 
tary transformations in a time-optimal way using avail- 
able physical resources is also a problem in quantum 
information processing. It has received significant at- 
tention, and time-optimal control of two coupled qubits 
0, S S 0, S, IE 0, MM is now wel1 understood. Recently, 
this problem has also been studied in the context of lin- 
early coupled three-qubit topologies [IE EL HI] , where 
significant savings in the implementation time of trilin- 
ear Hamiltonians and logic gates between indirectly cou- 
pled qubits were demonstrated over conventional meth- 
ods. Many of these ideas have found applications in 
efficient ways to propagate coherences along Ising spin 
chains [13, [3] • However, the complexity of the general 
problem of time-optimal control of multiple qubit topolo- 
gies is only beginning to be appreciated. 

In this article, we study the problem of finding the 
shortest pulse sequences for creating desired coherences 
in a linear chain of spins (of spin 1/2) weakly coupled 
with unequal Ising couplings. In particular, we start by 
considering the case of three linearly coupled spins and 



generalize our methods to linear spin chains. This study 
has immediate applications to multi-dimensional nuclear 
magnetic resonance (NMR) spectroscopy [la ]. In multi- 
dimensional NMR experiments [la] , starting from the 
thermal state of the spins, a multiple quantum coherence 
between spins is synthesized, which helps to correlate the 
frequencies of various spins. Efficient pulse sequences for 
creating such coherences help to improve the sensitivity 
of the experiments. 

One approach for solving these problems reduces the 
efficient synthesis of multiple spin order to the geometri- 
cal question of finding shortest paths on the sphere under 
a special metric 0, H, UU, [l6j] ■ The metric enforces the 
constraints on the quantum dynamics that arise because 
only limited Hamiltonians can be realized [I?], EE Ell • 
Such analogies between optimization problems related to 
steering dynamical systems with constraints and geom- 
etry have been well explored in areas of control theory 
[IE IH| and sub-Riemannian geometry [22j] . In this pa- 
per, we study in detail the metric and the geodesies that 
arise from the problem of efficient synthesis of multiple 
spin coherences between qubits (spin 1/2) that are indi- 
rectly coupled via unequal couplings to a third qubit. 



II. OPTIMAL CONTROL OF A LINEAR 
THREE-SPIN SYSTEM WITH UNEQUAL 
COUPLINGS 

We consider a linear chain of three spins placed in a 
static external magnetic field in the z direction with Ising 
type couplings between next neighbors [IE H3|- ^ n a 
suitably chosen (multiple) rotating frame which rotates 
with each spin at its resonant frequency the Hamiltonian 
that governs the free evolution of the spin system is given 
by the coupling Hamiltonian 



'Electronic address: haidong@mit.edu 
^Electronic address: zeier@cecs.harvard.edu 
tElectronic address: navin@hrl.harvard.edu 



H c — 2Ji2ll z l2z + ^Jl^llzHz- 



We use the notation lg v = l a . , where a,j = v for 
j = £ and dj = otherwise (see [HI). The matrices 



2 



/2, are 
2- 



/«:=(S5)/2 s /„:=(2 o*)A and 
the Pauli spin matrices and Iq ■— (J?) is the 2 x 
dimensional identity matrix. 

If the Larmor frequencies of the spins are well sepa- 
rated, each spin can be selectively excited by an appro- 
priate choice of the amplitude and phase of the radio- 
frequency (RF) field at its resonance frequency. The 
goal of the pulse designer is to derive explicit controls 
for the variables comprising of the frequency, amplitude 
and phase of the external RF field to effect a net unitary 
evolution U(t) most efficiently. 

We begin with the problem of finding the shortest pulse 
sequence that transform the initial polarization I\ z on the 
first spin to a multiple quantum coherence, i.e., 

hx — ► ^hzhzhz- 

Example 1. The conventional strategy for achieving 
this transfer is given by the following stages 



In the first stage of the transfer, operator I\ x evolves to 
Zliyhz under the natural coupling Hamiltonian H c in 
T\ = 7r/(2Ji2) units of time. This operator is then ro- 
tated to 2Ii y l2x by applying a hard (ir/2) y pulse on the 
second spin, which evolves to ^I\yI 2y I 3 z under the nat- 
ural coupling, in r 2 = tt/(2Jxs) units of time. Finally, 
hard (tt/2) x pulses on first and second spin prepare the 
desired final state. The total evolution time is then sim- 
ply n + t 2 = tt/(2J 12 ) + 7r/(2J ls ). 

We now study time-optimal designs for achieving this 
(and more general) transfers. To simplify notation, we 
introduce the following symbols (O) := Tr (Op), where 
Tr denotes the trace, for the expectation values of oper- 
ators O. Let xi = (h x ), x 2 = (2Ii y I 2z ), x 3 = (2I ly I 2x ), 
and X4 = (^h y I 2 yhz) and X = (xi, x 2 , x 3 , £ 4 ) T . By ex- 
pressing the time t, in units of 1/ J\ 2 , the evolution of the 
vector X is given by 



dX 
~dt 



'0-10 

1 — u 

u -k 

K k 



X = B(u,k)X, 



(1) 



where k — J 23 / J\ 2 and u = u(t) is the control param- 
eter representing the amplitude of the y-pulse on the 
second spin. Now, the problem of optimal transfer is 
to find the optimal u(t) for steering the system from 
(1,0,0, 0) T to (0,0, 0, 1) T in minimal time. We consider 
also the more general case of transfering the system from 
(cosa,sina,0,0) T to (0, 0, cos/3, sin/3) T , where we con- 
sider different a,(3 £ [0,7r/2]. 

Example 1 (continued). In this picture, the conven- 
tional method to transfer (1,0, 0,0) T to (0, 0,0, 1) T is 
described by setting u(t) = for t\ units of time, trans- 
fering (1, 0, 0, 0) T to (0, 1, 0, 0) T . Using the control u, we 
rotate (0, 1, 0, 0) T to (0, 0, 1, 0) T in arbitrary small time, 
as the control can be performed much faster as compared 



FIG. 1: Auxiliary variables r-2 and 0. 



to the evolution of couplings. Then, we set u(t) = and 
evolve t 2 units of time under the coupling Hamiltonian, 
transfering (0, 0, 1, 0) T to (0, 0, 0, 1) T . The total time for 
the transfer is n + t 2 . 

We now show how significantly shorter transfer times 
are achievable, if we relax the constraint that the selec- 
tive rotations on second spin are only hard pulses. We 
show that if we let the selective operations be carried out 
by soft shaped pulses, along with evolution of the cou- 
pling Hamiltonian, then we can achieve shorter transfer 
times. The pulse shapes can be numerically computed 
by formulating the problem as a derivation of geodesies 
on a sphere under a special metric as detailed below. 

We first make a change of variables (see Fig. [1]) . Let 



T\ = xi,r 2 



x 2 +^3' r 3 



X4 , and tan = — . 

x 2 



Using u(t), we can control the angle 9, so we can think 
of 6 as a control variable. Expressing the time in units of 
1/^12, the evolution of the system w.r.t. the coordinates 
Ti is given by 







I \r 2 ) = Icos^t) 

\rz) \ fcsin0(t) 



cos0(i) \ fn 
-ksm9(t) r 2 
/ U, 



_ , (2) 

The problem of transferring the system in Eq. ([T]) from 

(cosa,sina,0,0) T to (0, 0, cos/3, sin/3) T , reduces to find- 
ing 9(t), for steering the system from (cos a, sin a, 0) T 
to (0, cos (3, sin (3) T in minimal time in Eq. (J5J. Wc 
show that this is equivalent to finding the correspond- 
ing geodesic on the sphere, under the metric 



9 



k 2 dr\ 



drl 



k 2 r\ 



(3) 



By substituting for sin#(i) and cos9(t) from Eq. ^j, 
the transfer time t = J Q r -\/[sin 9(t)] 2 + [cos6'(i)] 2 dt re- 
duces to 



1 

T= k 



^[k 2 (h) 2 + (h) 2 }/r 2 2 dt=\ 



L dt. 



Thus, minimizing r, amounts to computing the geodesic 
under the metric g. The Euler-Lagrange equations for 
the geodesic take the form 



d_ t dL_ 
Jt V&h 



dL 
dri 



dt \dr 3 J V dr 3 



3 



Note, r\ = 1 — r\ — r\ and along the geodesic, L — k is 
constant. We get, 

I** and *(t*) =L *!* t ( 4 ) 
dt \ r\ ) r\ dt \rl) r 2 ' v 1 



which implies that 

d ( k 2 f\T 3 — r$r\ \ 



dt 



a,* - 1) 



nr 3 



(5) 



Let / = (fc 2 fir 3 — r 3 r\)jr\. From Eq. Q, we get 

dt \ r 2 J r 2 

Substitute with f\jr% = — cos 9(t) and f 3 /r 2 = fcsin0(i), 
and get 



[-fc 2 cos0(i)] = -/fcsin0(i), 



so = —f/k. Differentiating again, we obtain from 
Eq. ©, that 

= -14- f = - P ~ lh i 3 = (k 2 - f) cos0sin0. 
kdt J k r\ y ' 

Using a = k 2 — 1 , we rewrite this as, 



9 = ~sin20(t). 



(6) 



The solution to Eq. ([S]), can be given in terms of an 
elliptic integral. Note that by multiplying both sides of 
Eq. with 9, we get the equation of an ellipse in terms 
of the coordinates (0, cos[0]) T as 

c = 9 2 (t) + a cos 2 9 (t), (7) 

where c is a constant. Thus, we obtain that, 9 — 
±yjc — a cos 2 9(t) and 

r»V da 



(o) yc — a cos 2 a 



= ±t. 



(8) 



The left hand side of Eq. © is an elliptic integral [HH^. 
Equation ©, can be integrated if 9(0) and ^(0) are both 
known explicitly. In transfers considered subsequently, 
only 9(0) is known and therefore one has to numerically 
search over the possible values of 9(0), such that the 
resulting trajectory X(t) achieves the desired transfer. 
Note, guessing an initial value of 9(0) is same as search- 
ing for the correct value of c in Eq. ([7]) . 



III. COMPUTATION OF OPTIMAL CONTROL 
LAWS 

Now, consider the problem of steering the system of 
Eq. ©, from the initial state (r 1 ,r 2 ,r 3 ) T = (1,0, 0) T to 
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FIG. 2: For different k, we plot (a) the minimal time T in units 
of I/J12 and (b) the corresponding ratio T/[n/2 + 7r/(2fe)]. 



the target state (0, cos/3, sin/3) T . For j3 — n/2, we get 
the target state (0, 0, 1) T . Recall that, 



k6 



f = —k 2 cos 9— — ksva.9—. 



(9) 



At t = 0, n(Q) = 1 and r 2 (0) = r 3 (0) = 0. For / to be 
finite, at t = 0, we should have sin6>(0) = 0. Therefore, 
in Eq. ([6]), the initial condition is 9(0) = 0. To solve 
Eq. ©, we only need to know the initial value of 9(0), 
which is the same as knowing the constant c in Eq. I|7|). 
Given 9(0) = and the value of 0(0), we can numeri- 
cally solve first Eq. ^ and then Eq. ©. Consequently, 
we can determine for each value of 0(0), the smallest 
time s such that r\ (s) = by using a one-dimensional 
search, e.g., using the bisection method (see pp. 46-51 of 
Ref. [13 )• Thus, we can determine the value of 0(0) such 
that r 2 (s) — cos(/3), again by using a one-dimensional 
search. In summary, we have reduced the original opti- 
mzation problem to (combined) one-dimensional search 
problems. 

For different values of k, we plot the minimal time 
T = T(k) to transfer (n,r 2 ,r 3 ) T = (1, 0, 0) T to (0, 0, if 
in Fig.[2ja) in units of 1/J\ 2 . For example, when k = 1, 
i.e., J12 = J23 = J, it takes 2.72/ J amount of time, 
which is about 86.6% of tv/J, the time needed using the 
conventional method. Figure[2jb), shows the ratio of the 
the minimal time T(k) to the time ir/2 + ir/ (2k) obtained 
using the conventional strategy. We define this ratio rj(k). 
In both the plots k > 1 is considered, as the case for 
k < 1 can be derived from this. Figure [3l shows the 
optimal control u(t) in Eq. (fTJ) for different values of k. 

Observe that T(l/k) = kT(k). Let u(t,k) be the 
time optimal control for steering the system of Eq. (fl]) 
from (1, 0, 0, 0) T to (0, 0, 0, 1) T . Then the control v(t) = 
u(T — t, k), will steer the same system from (0, 0, 0, 1) T 
to (1, 0, 0, 0) T in the same time, which is also minimal for 
this transfer. Let Y = (x±, x 3 , x 2 , xi) T and consider the 
control v(t) — u(T — t,k). Then, we have 



dY 
~d~r 



kB(l/k,v)Y. 



But, we have just remarked that the minimal time to 
steer Y from (1,0,0,0) T to (0,0,0,1) T is T(k), and it 
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k = 1/10./3 = tt/2 




FIG. 3: The figure shows the shape of the optimal control 
u(t) in Eq. ([lj for k = 1, 2, 10. The time is in units of I/J12. 
The control amplitude is in units of J12 . Rightmost datapoint 
of each control u(t) is not displayed for reasons of numerical 
accuracy. 



follows that T(l/k)/k = T(k). It also follows that the 
optimal control u(t, 1/fc) = u[T(k) — t/k, k), where t £ 
[0, kT(k)]. Note t](k) = rj(l/k). 

Remark 1. In Fig. [U we present the geodesies of the 
metric of Eq. © for = tt/2 and k £ {1/10,1,10}, 
where (3 = tan" 1 [r 3 (T)/r 2 (T)], i.e., n(T) = r 2 (T) = 
and r 3 (T) — 1. Observe that the geodesies for the case 
k > 1, bend more towards the point (ri, r2, ra) = (0, 1, 0), 
before approaching the final point (0, 0, 1). In an intuitive 
way to understand this, consider first the limit k 1. 
In this limit, the minimum time to steer Eq. © from 
(1, 0, 0) to (0, 0, 1) is essentially the time required to steer 
the system (of Eq. |2])) to the equator. This is achieved 
fastest by moving from (1,0,0), directly to (0,1,0) and 
keeping 0(0) — 0. As k is decreased to 1, the geodesies 
gradually move away from (0, 1,0), implying 0(0) for k > 
1 is smaller than for the case k = 1. In the special case 
of k = 1, Eq. © reduces to = C, a constant. This case 
was been studied in detail in pjl, [l4| ■ Therefore, in the 
numerical computation of optimal control, the search for 
true 0(0) can be restricted to the interval [0, C]. 

In Fig. [5l we plot the minimal time T = T(k,f3) 
to transfer (ri,r 2 ,r 3 ) T = (1,0, 0) T to (0, cos/3, sm(3) T , 
for different values of k and (3 = tan _1 [r3(T)/r2(T)]. 
Figure O shows the geodesies of the metric of Eq. © 
for k = 2, for two different values (3 £ {7r/4, 7r/2} of the 
terminal point (0, cos /3, sin f3) T . 

If we transfer the system from (xi, x 2 , x 3 , X4) T = 
(cosa,sina,0,0) T to (0, 0, cos/3, sm(3) T , where a > 0, 
we do not know the value of 9(0). But from Eq. ([9]), 
we have the relationship 0(0) = sin[0(O)] cot(a) between 
0(0) and 0(0). Thus, we can apply similar search meth- 
ods for the right initial conditions as in the case of 



FIG. 4: Geodesies to transfer from (n, r^, r- A ) T = (1, 0, 0) T to 
(0,0, 1) T , where k = 1/10, k = 1, and k = 10. 




FIG. 5: For different k and (3 = tan" 1 [r 3 (T)/r 2 (T)], we plot 
the minimal time T in units of 1/ J12 



a = 0. For k = 2 and different values of a and (3, we 
plot in Fig. [7] the minimal time T = T(a, f3) to transfer 
(ri,r 2 ,r 3 ) T = (cos a, sin a, 0) T to (0, cos/3, sin/3) T . 



IV. EFFICIENT CREATION OF COHERENCES 
IN ISING SPIN CHAINS WITH UNEQUAL 
COUPLINGS 

We now consider a generalization of the problem 
treated above. We consider a linear chain of n spins, 
placed in a static external magnetic field in the z- 
direction, with unequal Ising type couplings between next 
neighbors 0, [24| . In a suitably chosen (multiple) rotat- 



■5 



k = 2,p = n/2 



l\ x to a multiple spin order represented by the opera- 



(1,0,0) 




FIG. 6: Geodesies to transfer from (n, r 2 , r 3 ) T = (1, 0, 0) T to 
(0,cos/3,sin/3) T , where k = 2 and (3 € {7r/4,7r/2}. 
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(rad) 



a (rad) 



FIG. 7: For fc = 2 and different values of a and /3, we plot 
the minimal time T in units of I/J12 



ing frame, which rotates with each spin at its resonant 
frequency, the free evolution of the spin system is given 
by the coupling Hamiltonian, 



H, 



i+ihzl(i+i) z - 



As before, we assume that the Larmor frequency of the 
spins are well separated compared to the coupling con- 
stant Jjjj+i, so that we can selectively rotate each spin at 
a rate much faster than the evolution of the couplings. 

We first consider the problem of synthesizing a uni- 
tary transformation U, which efficiently transfers a sin- 
gle spin coherence represented by the initial operator 



tor 2 n - i nr=i^ 



The conventional strategy for achieving this transfer is 
achieved through the following stages 



hx — » ^hyhz — > ^h y hx — » ^hyhyhz 
I n ~ 1 

^hyhyhx 2" 1 ( hy)Inz 

1=1 

where each evolution represents an appropriate evolution 
by rotation angle tt/2. The final state 2 n_1 (n£=i h y )Inz, 
is locally equivalent to 2™ _1 fl/Li The whole transfer 
involves n — 1 evolution steps under the natural Hamilto- 
nian, where the first step takes only 7r/(2J/ ; / + i) units of 
time, resulting in a total time of Y^i=i n /(2Ji,l+i)- We 
now formulate the problem of this transfer as a problem 
of optimal control and derive time efficient strategies for 
achieving this transfer. 

To simplify notation, we introduce the following sym- 
bols for the expectation values of operators that play a 
key part in the transfer: 

xi = (h x ),x 2 = (2I ly I 2z ), 
X3 = {2hyhx),x A = (4I ly I 2 yhz), ■ ■ ■ , 
X2n-3 = (2™ 1 hyhylzy ' ' ■I{n-l)x)i 
and X2n-2 = (2 n hyl2yl3y ' ' ' I(n-l)ylnz) ■ 

Let X — (xi,X2,xz, . . . ,X2n-2) T ■ Expressing the time 
in units of I/J12, the evolution of the system is given by 



dX_ 
~dt 



/0 -fci 

fcl —Ml 

Ml 











-k 2 

k 2 -U2 

m 2 -k 3 

fc 3 



x, (10) 



where ui are the control parameters representing the am- 
plitude of the y pulse on spin I + 1 and fc/ = J/,/+i/Ji2- 
The problem now is to find the optimal ui(t), steering 
the system from (1, 0, 0, • • • , 0) T to (0, 0, • • • , 0, 1) T in the 
minimal time. 

We divide the optimal transfer problem into mul- 
tiple steps. Let I = l,...,n — 1. Consider 
the step that steers (x 2 /_i, x 2 i, X21+1, x 2 i+2) T from the 
initial state (cos/3/, sin/3/, 0, 0) T to the target state 
(0, 0, cos/3/ + i, sin/3/ + i) T , by optimal choice of u/. 



(11) 



We denote the minimal time for this transfer as 
T(/3i, ki, fc/+i), as it depends on the parameters, 
{PuPi+ukuh+i)- Note that 

fc/T(/3/,/3/+i,fc/,fc/+i) = T(/3/,/3 ;+ i, l,fc/ + i/fc/). 





fX 2 l-l\ 
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FIG. 8: The figure shows the plot of function J (7) in Eq. (fl3jl . 
as a function of 7. 

We can then break the original transfer in Eq. (|TD|) . 
into sub-transfers denoted by 

fc-tfo > 01 ► 

where both 0\ = and /3 n _i = tt/2. The control in 
Eq. (jTTJ) should be so chosen that the time for transfer 
Pi to fli+x is minimized. Furthermore, the choice of in- 
termediate Pi should be made to minimize the total time 
of transfer 

n-2 

2=1 

The notation for Ji(0) denotes that the index I in the 
above summation begins with I — 1 with initial angle 
pi = 0. Now suppose, we have computed J;+i(£) for 
all £ 6 [0, 7r/2]. Then we can compute Ji(Pi) from the 
equation below as 

J i (A)=nain{J m (j7) + T(A ) ry I fe i ,fe i+ i)}. (12) 
v 

We can solve Eq. (Tr2"|) iteratively starting from 2 = n— 2 
and proceeding backwards and remembering that 

Jn-2{Pn-2) = T(/3 n _ 2 , 7t/2, fc„_ 2 , fc„-l). 

Note T(/3;, 77, fcz, fc;+i) in Eq. (fT2|) can be explicitly com- 
puted, as we showed in the previous section. Therefore 



Eq. (JTSJ) helps us to compute the optimal 77 and the cor- 
responding optimal control it;. 

Now, the problem is reduced to a simple dynamic pro- 
gramming problem, which can be solved efficiently. 

Example 2. We consider as an example, an Ising spin 
chain with four spins with Ji2/(27r) = 91 Hz, J 2 3/(27r) = 
15 Hz and Jm/(2tt) — 55 Hz. This system represents 
the popular [28[ HNCACO experiment in multidimen- 
sional NMR, where first spin is the proton, the second 
one represent 15 N, the third and fourth one are 13 C. For 
the transfer in Eq. (|10[) . the dynamic programming equa- 
tions for the minimum time is in units of I/J23, can be 
written as min 7 J (7), where 

J( 7 ) = {T(0, 7, fci, 1) + T( 7) |, 1, fc 2 )} (13) 
where fci = J12/J23 and fc 2 = ^34/^23- Figure [5] shows 
the plot of J( 7 ) and the minimum is achieved at 7op t ~ 
.1937T and the corresponding minimum value of J( 7op t) 
is » 2.01. The conventional transfer strategy will take 
7t(l + 1/fci + l/fc 2 )/2 w 2.26, units of time, which is 
w 12.2% longer than the proposed efficient methodology. 



V. CONCLUSION 

In this article, we studied the problem of efficient cre- 
ation of multiple spin order in an Ising spin chain with 
unequal couplings. We first analyzed in detail the system 
of three linearly coupled spins. We showed that the time 
optimal pulse sequences for creating multiple spin order 
in this system can be obtained by computing geodesies on 
a sphere under a special metric and that the solution to 
the resulting Euler Lagrange equation is related to the 
solution of a nonlinear oscillator equation 9 = AsmO. 
We showed that the minimum times and optimal control 
laws for the studied problem can be explicitly computed 
and provide significant gains over conventional methods. 
The methods developed in the three spin case were then 
exploited to find efficient strategies for manipulating the 
dynamics of Ising spin chains with unequal couplings. It 
is expected that these methods will find immediate ap- 
plications in coherent spectroscopy and quantum infor- 
mation processing. 
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